function xOut=plus(x1,x2)
  n1 = size(x1);
  n2 = size(x2);
  if(prod(n1)==1 & prod(n2)>1)
      x1 = repmat(x1,n2(1),n2(2));
  end
   if(prod(n2)==1 & prod(n1)>1)
      x2 = repmat(x2,n1(1),n1(2));
  end
  x = getval(x1);
  y = getval(x2);
  dfdx1 = ones(size(x1));
  dfdx2 = ones(size(x2));
 
  if(notisa(x1,'deriv1f'))
    np = nindep(x2);
  else
    np = nindep(x1);
  end
  if(notisa(x1,'deriv1f'))
    xOut.v = plus(x1,x2.v);
    xOut.d = elelmult_lastindx(dfdx2,x2.d);
  elseif(notisa(x2,'deriv1f'))
    xOut.v = plus(x1.v,x2);
    xOut.d = elelmult_lastindx(dfdx1,x1.d);
  else  % both deriv1f
    n1d = size(x1.d);
    np = n1d(end);
    xOut.v = plus(x1.v,x2.v);
    xOut.d = elelmult_lastindx(dfdx1,x1.d) + elelmult_lastindx(dfdx2,x2.d);
  end
  xOut=class(xOut,'deriv1f');
